{
 "metadata": {
  "name": ""
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import numpy as np\n",
      "from Factor import *\n",
      "from FactorOperations import *\n",
      "from PGMcommon import *\n",
      "import profile\n",
      "from CliqueTree import *\n",
      "from CliqueTreeOperations import *\n",
      "import networkx as nx\n",
      "import matplotlib.pyplot as plt\n",
      "\n",
      "\"\"\" Lets try and construct a basic Bayesian network based on the Student example given\n",
      "in Koller and Friedman on pg. 53 (Figure 3.4)\n",
      "\n",
      "Made the Student network in samiam software <http://reasoning.cs.ucla.edu/samiam/help/>\n",
      "saved the file to /Users/amit/software/pgmPy/Student.net.\n",
      "\n",
      "Read in matlab like so: studentBN=ConvertNetwork('/Users/amit/software/pgmPy/Student.net')\n",
      "Then was able to create the Factor objects like so below:\n",
      "\"\"\"\n",
      "\n",
      "D=Factor( [1], [2], [.6,.4],'Difficulty')\n",
      "I=Factor( [5], [2], [.7,.3], 'Intelligence' )\n",
      "G=Factor( [2, 5,1], [ 3, 2, 2], [.3,.4,.3,.9,.08,.02,.05,.25,.7,.5,.3,.2], 'Grade', )\n",
      "S=Factor( [3,5], [2,2], [.95,.05,.2,.8], 'SAT')\n",
      "L=Factor( [4,2], [2,3], [.1,.9,.4,.6,.99,.01], 'Letter' )\n",
      "\n",
      " \n",
      "#print np.prod( G.getCard() )\n",
      "#print IndexToAssignment( np.arange(np.prod(G.getCard())), G.getCard() )\n",
      "\n",
      "INPUTS=[D,G,S,L,I] # list of Factors that comprise the Bayesian Network\n",
      "                   # see StudentBN.jpg\n",
      "jointFactor = ComputeJointDistribution(INPUTS)\n",
      "#print jointFactor\n",
      "\n",
      "\n",
      "#print ExpFactorNormalize( FactorSum( LogFactor(D), LogFactor(G) ) )\n",
      "#print FactorProduct(D, G)\n",
      "#print np.sum( jointFactor.getVal() )\n",
      "\n",
      "print jointFactor.getVar()\n",
      "print jointFactor.getCard()\n",
      "Ia=IndexToAssignment ( np.arange(np.prod(jointFactor.getCard()) ), jointFactor.getCard() ) -1\n",
      "print Ia\n",
      "print \"===\"\n",
      "#print Ia[:,4]\n",
      "print Ia[34]-1\n",
      "#print jointFactor.getVal()[34]\n",
      "#print GetValueOfAssignment(jointFactor, np.array([[0,2,1,0,1]]))\n",
      "#print len(jointFactor.getVal())\n",
      "#IndexToAssignment(1:np.prod(P.cliqueList(1).card), P.cliqueList(1).card)\n",
      "\n",
      "#elimination=[3,5,1,4,2]\n",
      "\n",
      "#print MaxProductVE(elimination, INPUTS)\n",
      "#print ComputeMarginal( [4], INPUTS, [])\n",
      "\n",
      "\n",
      "#clique_beliefs=ComputeExactMarginalsBP ( INPUTS)\n",
      "#for c in clique_beliefs:\n",
      "#    print c\n",
      "#    print"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[1 2 3 4 5]\n",
        "[ 2.  3.  2.  2.  2.]\n",
        "[[ 0.  0.  0.  0.  0.]\n",
        " [ 1.  0.  0.  0.  0.]\n",
        " [ 0.  1.  0.  0.  0.]\n",
        " [ 1.  1.  0.  0.  0.]\n",
        " [ 0.  2.  0.  0.  0.]\n",
        " [ 1.  2.  0.  0.  0.]\n",
        " [ 0.  0.  1.  0.  0.]\n",
        " [ 1.  0.  1.  0.  0.]\n",
        " [ 0.  1.  1.  0.  0.]\n",
        " [ 1.  1.  1.  0.  0.]\n",
        " [ 0.  2.  1.  0.  0.]\n",
        " [ 1.  2.  1.  0.  0.]\n",
        " [ 0.  0.  0.  1.  0.]\n",
        " [ 1.  0.  0.  1.  0.]\n",
        " [ 0.  1.  0.  1.  0.]\n",
        " [ 1.  1.  0.  1.  0.]\n",
        " [ 0.  2.  0.  1.  0.]\n",
        " [ 1.  2.  0.  1.  0.]\n",
        " [ 0.  0.  1.  1.  0.]\n",
        " [ 1.  0.  1.  1.  0.]\n",
        " [ 0.  1.  1.  1.  0.]\n",
        " [ 1.  1.  1.  1.  0.]\n",
        " [ 0.  2.  1.  1.  0.]\n",
        " [ 1.  2.  1.  1.  0.]\n",
        " [ 0.  0.  0.  0.  1.]\n",
        " [ 1.  0.  0.  0.  1.]\n",
        " [ 0.  1.  0.  0.  1.]\n",
        " [ 1.  1.  0.  0.  1.]\n",
        " [ 0.  2.  0.  0.  1.]\n",
        " [ 1.  2.  0.  0.  1.]\n",
        " [ 0.  0.  1.  0.  1.]\n",
        " [ 1.  0.  1.  0.  1.]\n",
        " [ 0.  1.  1.  0.  1.]\n",
        " [ 1.  1.  1.  0.  1.]\n",
        " [ 0.  2.  1.  0.  1.]\n",
        " [ 1.  2.  1.  0.  1.]\n",
        " [ 0.  0.  0.  1.  1.]\n",
        " [ 1.  0.  0.  1.  1.]\n",
        " [ 0.  1.  0.  1.  1.]\n",
        " [ 1.  1.  0.  1.  1.]\n",
        " [ 0.  2.  0.  1.  1.]\n",
        " [ 1.  2.  0.  1.  1.]\n",
        " [ 0.  0.  1.  1.  1.]\n",
        " [ 1.  0.  1.  1.  1.]\n",
        " [ 0.  1.  1.  1.  1.]\n",
        " [ 1.  1.  1.  1.  1.]\n",
        " [ 0.  2.  1.  1.  1.]\n",
        " [ 1.  2.  1.  1.  1.]]\n",
        "===\n",
        "[-1.  1.  0. -1.  0.]\n"
       ]
      }
     ],
     "prompt_number": 122
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "P=CreatePrunedInitCtree(INPUTS, [] )\n",
      "(P, MESSAGES) = CliqueTreeCalibrate(P, isMax=1)\n",
      "print P.getEdges()\n",
      "#cliqueList=P.getNodeList()\n",
      "#for i in range( len(cliqueList)):\n",
      "#    print cliqueList[i]\n",
      "#    print\n",
      "#print \"===\"\n",
      "\n",
      "\n",
      "#print FactorMaxMarginalization( cliqueList[0], [3] )\n",
      "#print FactorMaxMarginalization( cliqueList[2], [1,2] )\n",
      "#G=nx.from_numpy_matrix( P.getEdges() )\n",
      "#nx.draw_shell(G)\n",
      "#plt.show()\n",
      "#print \"====\"\n",
      "#cliqueNodes=P.getNodeList()\n",
      "print MESSAGES[0,2].getVal()\n",
      "print np.exp( MESSAGES[0,2].getVal() )\n",
      "print \"====\"\n",
      "#print lognormalize( MESSAGES[2,0].getVal() )\n",
      "print MESSAGES[2,0].getVal()\n",
      "print np.exp( MESSAGES[2,0].getVal() )\n",
      "#print \"====\"\n",
      "\n",
      "jointDistribution=ComputeJointDistributionFromCalibratedCliqueTree(P, MESSAGES, 1)\n",
      "print jointDistribution\n",
      "print np.sum( jointDistribution.getVal()/np.sum(jointDistribution.getVal()) )\n",
      "#pSNP=jointDistribution.getVal()[0]\n",
      "#print 1-pSNP"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[[ 0.  0.  1.]\n",
        " [ 0.  0.  1.]\n",
        " [ 1.  1.  0.]]\n",
        "[-0.40796824 -1.42711636]\n",
        "[ 0.665  0.24 ]\n",
        "====\n",
        "[-1.28301601 -0.72154666]\n",
        "[ 0.2772  0.486 ]\n",
        "name: None\n",
        "var: 1 2 3 4 5\n",
        "card: 2.0 3.0 2.0 2.0 2.0\n",
        "val: 0.113847191769 0.0126496879744 0.607185022769 0.252993759487 1.12708719852 1.75324675325 0.00599195746154 0.000665773051282 0.0319571064615 0.0133154610256 0.0593203788693 0.0922761449077 1.02462472592 0.113847191769 0.910777534154 0.379490639231 0.0113847191769 0.0177095631641 0.0539276171539 0.00599195746154 0.0479356596923 0.0199731915385 0.000599195746154 0.000932082271795 0.0308157812308 0.0114132523077 0.0109567222154 0.0273918055385 0.00677947187077 0.0451964791385 0.123263124923 0.0456530092308 0.0438268888616 0.109567222154 0.0271178874831 0.180785916554 0.277342031077 0.102719270769 0.0164350833231 0.0410877083077 6.84795138462e-05 0.000456530092308 1.10936812431 0.410877083077 0.0657403332923 0.164350833231 0.000273918055385 0.00182612036923\n",
        "1.0\n"
       ]
      },
      {
       "output_type": "stream",
       "stream": "stderr",
       "text": [
        "A factor is empty!\n",
        "A factor is empty!\n",
        "A factor is empty!\n",
        "A factor is empty!\n",
        "A factor is empty!\n",
        "FactorMarginalization:Error, resultant factor has empty scope...\n"
       ]
      }
     ],
     "prompt_number": 123
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "\n",
      "\"\"\" factor division example Koller and Friedman page 365 \"\"\"\n",
      "N=Factor( [1,2], [3,2], [.5,0,.3,.2,0,.45],'Numerator')\n",
      "D=Factor([1], [3], [.8,0,.6], 'Denominator')\n",
      "print N\n",
      "print D\n",
      "print \"===\"\n",
      "print FactorDiv(N,D)\n"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "name: Numerator\n",
        "var: 1 2\n",
        "card: 3 2\n",
        "val: 0.5 0.0 0.3 0.2 0.0 0.45\n",
        "name: Denominator\n",
        "var: 1\n",
        "card: 3\n",
        "val: 0.8 0.0 0.6\n",
        "===\n",
        "name: None\n",
        "var: 1 2\n",
        "card: 3.0 2.0\n",
        "val: 0.625 0.0 0.5 0.25 0.0 0.75\n"
       ]
      }
     ],
     "prompt_number": 124
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "\"\"\" given a Factor object F, calculate its variable stride in value array of the factor \"\"\"\n",
      "\n",
      "def variableStride( f ):\n",
      "    strides=[]\n",
      "    iF=IndexToAssignment ( np.arange(np.prod(f.getCard()) ), f.getCard() )\n",
      "    variables=f.getVar()\n",
      "    cardinalties=f.getCard()\n",
      "    for i in range(0, len(variables) ):\n",
      "        assignment_slice=iF[:,i]\n",
      "        var_card=cardinalties[i]\n",
      "        #print variables[i], cardinalties[i]\n",
      "        curr_stride=iF[:,i]\n",
      "        if cardinalties[i] > 1:\n",
      "            stride=np.where(curr_stride==2)\n",
      "        else:\n",
      "            stride=0\n",
      "        #print variables[i], cardinalties[i], stride[0][0]\n",
      "        strides.append( stride[0][0] )\n",
      "        #print \n",
      "    \n",
      "    return strides"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 158
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "strides=variableStride( jointFactor )"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 159
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print strides"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[1, 2, 6, 12, 24]\n"
       ]
      }
     ],
     "prompt_number": 160
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "def IndexOfAssignment( f, strides, assignment):\n",
      "    \"\"\" given a Factor object f and the strides of each of the variables of f and \n",
      "    a list of assignments of those variables, return the position in the value array of that assignment\n",
      "    see page 358 box 10.a in Friedman and Koller text on how efficiently map a particular variable\n",
      "    assignment to the and index in the value array of a factor\n",
      "    \"\"\"\n",
      "    \n",
      "    idx=0\n",
      "    \n",
      "    cardinalities=f.getCard()\n",
      "    for (ass, stride) in zip( assignment, strides):\n",
      "        #print ass, stride\n",
      "        idx +=(ass *stride)\n",
      "\n",
      "    return idx\n",
      "    \n",
      "    "
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 161
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "idx_ass=IndexOfAssignment( jointFactor, strides, [0,2,1,0,1])\n",
      "print idx_ass\n",
      "\n",
      "idx_ass=IndexOfAssignment( jointFactor, strides, [0,0,0,0,0])\n",
      "print idx_ass"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "0 1\n",
        "2 2\n",
        "1 6\n",
        "0 12\n",
        "1 24\n",
        "34\n",
        "0 1\n",
        "0 2\n",
        "0 6\n",
        "0 12\n",
        "0 24\n",
        "0\n"
       ]
      }
     ],
     "prompt_number": 167
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print Ia[idx_ass]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[ 0.  0.  0.  0.  0.]\n"
       ]
      }
     ],
     "prompt_number": 168
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print jointFactor.getVal()[idx_ass]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "0.0028512\n"
       ]
      }
     ],
     "prompt_number": 164
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print jointFactor.getVal()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[ 0.01197    0.00133    0.06384    0.0266     0.118503   0.184338   0.00063\n",
        "  0.00007    0.00336    0.0014     0.006237   0.009702   0.10773    0.01197\n",
        "  0.09576    0.0399     0.001197   0.001862   0.00567    0.00063    0.00504\n",
        "  0.0021     0.000063   0.000098   0.00324    0.0012     0.001152   0.00288\n",
        "  0.0007128  0.004752   0.01296    0.0048     0.004608   0.01152    0.0028512\n",
        "  0.019008   0.02916    0.0108     0.001728   0.00432    0.0000072\n",
        "  0.000048   0.11664    0.0432     0.006912   0.01728    0.0000288\n",
        "  0.000192 ]\n"
       ]
      }
     ],
     "prompt_number": 135
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "val  = jointFactor.getVal().tolist()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 147
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print val"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[0.01197, 0.0013300000000000002, 0.06384, 0.026600000000000002, 0.11850299999999998, 0.18433799999999997, 0.0006299999999999999, 7.000000000000002e-05, 0.00336, 0.0014000000000000002, 0.0062369999999999995, 0.009701999999999999, 0.10772999999999998, 0.011970000000000001, 0.09575999999999998, 0.03989999999999999, 0.001197, 0.0018619999999999995, 0.00567, 0.0006300000000000001, 0.005039999999999999, 0.0021000000000000003, 6.299999999999999e-05, 9.799999999999998e-05, 0.0032400000000000007, 0.0012000000000000003, 0.001152, 0.00288, 0.0007128, 0.004752000000000001, 0.012960000000000003, 0.004800000000000001, 0.004608, 0.01152, 0.0028512, 0.019008000000000004, 0.029160000000000002, 0.010800000000000002, 0.0017280000000000002, 0.00432, 7.2000000000000005e-06, 4.800000000000001e-05, 0.11664000000000001, 0.04320000000000001, 0.006912000000000001, 0.01728, 2.8800000000000002e-05, 0.00019200000000000003]\n"
       ]
      }
     ],
     "prompt_number": 148
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print val.index(0.004608)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "32\n"
       ]
      }
     ],
     "prompt_number": 149
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "print Ia[32]\n"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "[ 0.  1.  1.  0.  1.]\n"
       ]
      }
     ],
     "prompt_number": 151
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [],
     "language": "python",
     "metadata": {},
     "outputs": []
    }
   ],
   "metadata": {}
  }
 ]
}